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DYNAMIC DATA ANALYSIS TECHNIQUES USED IN THE LANGLEY 
TIME SEMES ANALYSIS COMPUTER PROGRAM 


By Robert C. Ward 
Langley Research Center 

SUMMARY 

This paper is intended to be a reference guide to users of the Langley time series 
analysis program (designated R1299), The bases, for choosing available options, 
descriptions of the options, possible output quantities, equations of the various computa- 
tions, and useful tables and figures which might aid in the analysis of data are included. 

INTRODUCTION 

The Langley digital time series analysis (TSA) program is a general purpose pro- 
gram for the analysis of random, stationary time series, and the relationships between 
these series. The major output quantities are probability density estimators (amplitude 
domain), correlation estimators (time domain), and spectral estimators (frequency 
domain). A choice of standard (histograms for probability densities; correlations and 
Fourier transforms for spectra) or novel techniques (orthogonal functions for probability 
densities; recursive filtering for spectra) may be specified. 

Up to 20 different channels of data may be processed at one time. The input data 
may be conditioned by removal of its mean, corrected for a different bias in each channel, 
and subjected to a generalized preconditioning (shaping) filter with appropriate post- 
darkening of the spectra output. The preconditioning filter allows an arbitrary cascade 
of first-difference, high-pass, low-pass, band-pass, and band-reject filters. 

The actual processing of the conditioned data may include the computation of auto- 
covariances and cross-covariances; autocorrelations and cross-correlations; autospectral 
and cross -spectral densities; average frequency, equivalent bandwidth, and equivalent 
duration in each spectral band; coherences; interpowers; transfer functions; filter 
response functions; time to maximum cross-corr elation; spectral confidence bands; 
means; standard deviations; normalized third and fourth central moments; normality 
tests; histograms; and orthogonal coefficients for probability distribution. 

Since the TSA program is the major analytical tool for the analysis of dynamic data 
at Langley Research Center, a reference manual would be of considerable benefit to users 



or prospective users of this program. The purpose of this paper is to supply the user 
with this reference manual. (One should be aware of the different definitions placed on 
such terms as autocovariance mid autocorrelation which differ from one book to the next 
The defining equations should aiv/ays be checked before comparing results from two 
references.) 


Some of the information found in this paper was obtained from reference 1 , which 
is the documentation of the TSA parent program written under NASA Contract 
No. NAS 1-5632 and the information is included in this document for completeness . 


SYMBOLS 

A^^(f') filter response function for analysis filters 

Aj resultant nonrecursive weight 

a(^) cascaded nonrecursive weight for ith precondition filter 

■^ 0 (^ 0 ) quadrature filter nonrecursive weight 

•^0 (^ 0 ) in-phase filter nonre cursive weight 
aj coefficients for probability density function expansion 

*t nonrecursive weight for ith precondition filter In cascade 


Bi 

Bffe) 

b(i) 

m 

R 

pc 

b 


resultant recursive weight; bandwidth for frequencies in filtering technique 
arbitrary analysis 

quadrature filter recursive weight 

in-phase filter recursive weight 

cascaded recursive weight for ith precondition filter 

percentage confidence band 

equivalent bandwidth for channel x 

an endpoint of interval containing da±a samples 
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bi(fo) recursive weight to cascade for analysis filters 

recuxsive weight for ith precondition filter in cascade 

Cxy(f) real part of cross-spectral density between channels x and y 

(f) real part of normalized cross-spectral density between channels x and y 

xy 

D resultant number of recursive weights for precondition filter 

D(fo) decimation factor for analysis filter 

Di number of recursive weights for ith filter in cascade 

Dx(fo) equivalent duration for channel x 

di number of recursive weights for precondition filter 

dp precondition decimation factor 

Fj(z^^ Hermite or Laguerre orthonormal functions 

f frequency 

f' normalized frequency, f/fg 

Af frequency resolution for power spectra in correlation technique; filter 

bandwidth in filtering technique equispaced analysis 

fg filter bandwidth 

f' normalized filter bandwidth 

Jd 

f|, normalized precondition filter cutoff frequency 

fj maximum frequency of interest in analog data 

f| ith frequency to be analyzed in filtering technique 
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max 


maximum frequency to be analyzed in filtering technique 


% 

fo 

4 


f x(^o) 

Gx(f),Gy(f) 

G^(f) 

Gxy(f) 

H(f) 

Hf4(f’) 

Hx(f),Hy(f) 

I(x) 


T 

^xy 

K 

k 

M 


Nyquist frequency, fg/2 

filter center frequency 

normalized filter center frequency, fo/fs 

center frequency of filter preceding fg 

sampling frequency 

average frequency for channel x 

power spectral density for channels x and y 
normalized power spectral density for channel x 
amplitude of cross -spectral density between channels x and y 
instrumentation system amplitude response 
complex filter response 

amplitude frequency response for channels x and y from calibration tables 
integer part of x 

nth data point for channels x and y, respectively, filtered through 
quadrature filter with fo center frequency 

interpower between channels x and y 

fourth central moment (kurtosis); number of degrees of freedom 

index; number of filters per octave in filtering technique constant Q 
analysis 

number of filtered data points after decimation 
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m miBiber of lags in correlation technique 

N number of data samples to analyze; resultant number of nonrecursive weights 

Ni number of nonrecursive weights minus one for ith filter in cascade; number 

of points in ith histogram bin 

Nw number of filter warmup samples 

n^ number of nonrecursive weights minus one 

P number of expansion terms in probability expansion 

P(x) smooth probability density function 

P(z) smooth normalized probability density function 

Pj^ probability of being in ith histogram bin if Gaussian 

Qjjy(f) imaginary part of cross -spectral density between channels x and y 

(f) imaginary part of normalized cross -spectral density between channels x and y 

Rx(T),Ry(T) autocovariances of channels x and y, respectively 

Rxy(r) cross-covariance between channels x and y 

R^(fo),R^(fo) nth data point for channels x and y, respectively, filtered through 
in-phase filter with center frequency 

r harmonic number in power spectral density equation 

decibel spread for percentage confidence band 

T time length of data to be digitized 

T(f) precondition filter transfer function 

T^y(f) amplitude transfer function between channels x and y 



sampling interval (time increment) 


¥ 


number of analysis filter weights 


X 

X 

Xi 

Xi 


x,y 


yi 


zi 

a 

l3 

m 

^(fo) 

^x(f),/3y(f) 

y 

yxy(f) 


number of filters in cascade 

mean value 

bias for channel x 

input data for channel x with bias subtracted 

data sample with or without the mean and/or bias subtracted 

preconditioned data sample with or without the mean and/or bias subtracted 

data channels 

data sample for channel x 

data sample for channel y 

normalized data sample 

integration variable 

third central moment (skewness) 

instrumentation system phase response 

noise bandwidth of analysis filter pair 

phase response for channels x and y, respectively, from calibration tables 
intermediate variable; test for normality value 
coherence between channels x and y 

"lag window" (smoothing parameter) defined in power spectral density 
equation 
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(t) 


xy 


p^i'd 


T. 


xy 

p 

xy 




<!> (f) 

XyV / 


* (f) 

xy^ ^ 


%,a 

Subscripts: 

imag 

real 


max 


effective filter center frequency after decimation 
cross-correlation between channels x and y 

autocorrelation of channel x 

standard deviation 

lag time (displacement) 

time lag between channels x and y 

time lag to maximum cross-correlation between channels x and y 

analysis filter phase response 

analysis filter pair phase response difference 

phase of cross-spectral density between channels x and y 

phase of normalized cross -spectral density between channels x and 

chi-square value for K degrees of freedom and a type 1 error 

imaginary 

real 

maximum 


y 


DYNAMIC DATA REDUCTION 

The general procedure for reducing dynamic data in recorded analog form is found 
in figure 1. The solid lines represent required steps and the dashed lines represent 
optional steps and output. The digitizing process obtains a data sample every At sec- 
onds where At is determined by the sampling frequency as discussed later. The units 
of the digitized data are normalized counts. The quantity pass program converts the data 
units from normalized counts to the physical units as recorded. The time-history plot 
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and list program produces a listing of all or part of the input data and/or a "calcomp" 
plot of the data. The TSA output calcomp plot program produces calcomp plots of the 
TSA output functions. Documentation and information concerning these programs are 

available upon request from the Langley Analysis and Computation Division. 

If the data are already in digital form, the only step preceding the analysis per- 
formed by the TSA program is a format change. The data must be changed to the format 
accepted by the TSA program as shown in table I. If at all possible, the data should also 
be blocked as outlined in table II, although provisions are made to accept either blocked 
or unblocked data. An example of a FORTRAN WRITE statement for producing an 
unblocked data tape with six data channels is as follows: 

WRITE (1) (DATA(I),I=1, 16). 

DATA (3) through DATA (8) should have been equivalenced to a fixed-point array (dimen- 
sion 6) and the Hollerith engineering identification entered in this array. An example of 
a FORTRAN WRITE statement for producing a blocked data tape with six data channels 
is as follows: 

WRITE (1) (DATA(I),I=1, 500). 

DATA (1) through DATA (16), DATA (21) through DATA (36), etc., correspond to this 
example of the unblocked data array. 

Figure 2 shows a simplified flow of data through the three major sections of the 
TSA program. These sections are the probability density estimation, preconditioning of 
the data, and the spectral estimation. Each of these blocks as well as the computation 
method within each block is independent of the others. Thus, the user may choose not 
to compute the probability density function or choose either computation method desired 
without affecting his decision concerning computation of the spectral estimators or pre- 
conditioning of the data. 

SAMPLING FREQUENCY FOR DIGITIZING 

The sampling frequency fg at which the data is digitized is a very important fac- 
tor in the dynamic data reduction process. This frequency determines the Nyquist fre- 
quency fjq, which is the highest frequency that can be analyzed, by the following 
relationship; 
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Therefore, one could determine the necessary sampling frequency in theory by the equa- 
tion fg = 2f| where fj is the maximum frequency of interest contained in the data. 
With fg determined in this manner, there are only two points per cycle for determining 
the power at this maximum frequency of interest. Although two points per cycle of the 
frequency fj is the theoretical requirement, more points are recommended in practice 
for improved results. The suggested procedure at Langley Research Center is to select 
fg as follows: 

fg = 5% (2) 

Central data transcription at Langley Research Center automatically filters the 
data with an analog low pass filter before the data are digitized. This filter eliminates 
the problem of folding known as "aliasing." (See ref. 2, pp. 279-281.) 

In the remaining part of this paper, input data refer to the data after it has been 
digitized, converted to physical units, and recorded in the format shown in table I. 

TREATMENT OF THE INPUT DATA 

Several options are available for calibration and preconditioning of the input data. 
These options may be chosen independently of each other and of the methods chosen for 
the actual processing. The options are listed and the governing equations are given. 

Bias Calibration 

A bias level (zero level) value may be entered to subtract from each channel. The 
bias is subtracted before any computations are performed on the data; therefore, all out- 
put quantities reflect the data with the bias subtracted. If a bias is not entered for a 
particular channel, the bias to be subtracted for that channel equals zero. 

Xi = Xi - X (3) 


where 


Xi 

input data for channel x 

X 

bias for channel x 


Mean Removal 

By the proper selection of an input parameter, the mean may be removed from each 
channel before the probability and spectral estimators are computed. The mean is com- 
puted after the bias, if any, is subtracted from each data point . 
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If the mean is not removed, it will show up in the power spectral density function 


y-'i-P V-1.-1 rn vt vv% <-» r» 4 - -pv> ^ y-\vH 

tji lll'CcXlt liJ.£l^iAXvU.W*C CtC £aXZ:i. KJ Xi C:y£iac:Ai%-yy , 


tion theoretically is infinite and if there were any frequency components in the data near 
zero, one would not be able to determine them because they would be completely masked 

by the delta function. In actual practice, the delta function would just be a very sharp 
peak near zero frequency and would have a finite bandwidth. This condition is one of the 
reasons why the mean, in the general case, should be removed. 

Xi = Xi - X (4) 

where 


Xi input data for channel x with bias subtracted 

X mean value of channel x 


Preconditioning With a Shaping Filter 

Spectral computations at low power level are considerably distorted by the presence 
of high power in adjacent frequency bands. These distortions are reduced by the tech- 
nique of preconditioning the data by a shaping filter which tends to flatten the spectra or 
make it look more like white noise. (The effect of this filter may then be removed from 
the output by dividing the spectra by the filter transfer function. This procedure is 
called postdarkening and is discussed later under "Modification of Spectral Output.") The 
shaping filters may also be used to remove undesirable frequencies from the data before 
analysis. 

The shaping filter design allows any of the following different types of simple 
filters; first difference (an approximation to the first derivative filter), low pass, high 
pass, band pass, and band reject. Also, it allows an arbitrary series combination (cas- 
cading) of the preceding types. A filter made up of a series combination of simple fil- 
ters is called a compound filter. A preconditioning restriction is that all channels must 
be treated with the same filter. 

Figure 3 shows some typical transfer functions of the simple filters. Figure 4 
shows two examples of cascaded compound shaping filters. These transfer functions 
have an abscissa of normalized frequency f. Normalized frequency is the frequency 
divided by the sampling frequency. 

The filters are designed by use of a recursive technique. This technique is desir- 
able since it requires a smaller amount of computer storage and less time to perform 
the desired function. A recursive filter is a filter which has a feedback of the output 
samples into the filter. This filter is illustrated in figure 5 which shows both the pictoral 
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and mathematical representation. The two summations in the mathematical equation are 

the iionrecursive summation with the noiirecursive weights Ai a,nci the recursive sum- 

■ ^ J 

matioii with the recursive weights Bjn. The type of shaping filter designed is com- 
pletely determined by Aj, Bm, the number of nonrecursive weights N, and the number 
of recursive weights D. A note of caution in specifying these filters is that the recur- 
sive filter takes a time approximately equal to the reciprocal of the filter bandwidth to 
stabilize. If this time is more than 10 percent of the data span, the spectral estimators 
may be adversely affected. 


Preconditioning (if requested) is implemented by a shaping filter containing non- 
recursive weights Aj and recursive weights Bm. If a compound filter is requested, 
these weights are computed by cascading the weights of the simple filters. The ith simple 


filter has ni + 1 nonrecursive weights a^^\ 


.(i) 


and di recursive 


weights b 


(i) h(i) 


.(i) 


The various simple filter classes have the following synthesis; 


First difference: 


ni= 1 
di = 0 




A 


0.5 

- 0 . 


5 


Low pass (for cascading): 


Hi = 0 
di = 1 


A 


b 

a 


= exp 

(i)_i 
0 " 



( 5 ) 


( 6 ) 


Low pass (for simple filter); 


ni = l 
di = l 


.( 1 ) 


0.5 

.(1) 


1 + tan^rf^ 


sec 


(rfc) 


( 7 ) 
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High pass (for cascading); 


Hi = 0 


.(i) . 


= -exp^-27rfc^ 


^0 


4‘> 


( 8 ) 


High pass (for simple filter): 


n^ = l 




di = 1 


M) 


= O.sjl + sec^Trfg^ - tan^vrfc^ 


ag - -aj 

= 2a^j^^ - 1 


( 9 ) 




Band pass: 


nj^ = 0 


d^ = 2 


A 


b^^^ = -exp^-27rfg^ 
bj^^) = -2b^2^) cos(27rfo) 

+ (bW)^ + (b^^^ cos(27Tf’) - 2b^ cos(47rf') 


aW- 

ag - 


.(i) 


+ 2b^^^b^^^ cos^27rfQ^ 


1/2 


( 10 ) 


J 


This filter is a band-pass filter -with zero phase at the normalized center frequency 1 q 
and the maximum response slightly greater than 1, 



Band reject; 


di=2 


ii) 


-exp 


("2^b) 


= -2b2 ^ cos^27rfQ 


4‘> = -4*’ 

.<« = -b« 


a(i) - 1 - 
aQ - J- 


1 -f- b 


+ 2bj^^b2^^ cos(27rf') 


. (b«^' 


>1 y T ^^^2 y - 2b^^^ cos(27Ti;') - 2b^^^ cos 


(4^4) 


1/2 




( 11 ) 


The weights for the ith cascaded filter are built up from the (i - l)th filter by the fol- 
lowing recurrence relationships: 




k=0 


"k ^i-k 


(i = o, 1, 


Ni) 


( 12 ) 


g(i) = 


di 

I 

k=l 


b(i)B(i“l) 

°k ^m-k 


(m = 1, 2, . . Di) 


(13) 


where 


Ni = N._^ 4- n^ 


Di = D._j di 

with the initial conditions; 

a(o) = 1 

Aq 

= 0 

] 

■n(^) _ 


(i * 0) 


B„ 


" rn 


(m 0) 
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When X| is used to denote the ith data sample with or without the mean and/or bias 
subtracted, the application of the shaping filter is then: 


N 

D 

V — \ 


\ A^X. . + 

> B^x! 

™ 1 -m 

(14) 

J=0 

m=l 



1 

If the spectra are to be computed by the filtering technique and fmax the shaping 

5 

filter output is decimated by dp which is determined by the following equation: 

d =l/l— i_\ (15) 

where 

At time increment 

I(x) integer part of x 

fmax maximum frequency to analyze in filtering technique 

If preconditioning is requested, all spectral answers are computed by using X| and all 
probabilities are computed by using Xp In the remaining equations, x^ will be used 
for Xi or X^ depending on the treatment of the data. 

SPECTRAL METHODS 

There are two types of spectral computations available for use in the program. 
They are (1) the traditional correlation- Fourier transform techniques and (2) a novel 
filtering technique that is directly analogous to an analog filter bank. 

The correlation- Fourier transform method is to be preferred when any of the fol- 
lowing occurs: 

(1) Autocorrelation or cross -correlation functions are required for data reduction. 

(2) Cross-spectra are required on a small number of channels (small as compared 
with the number of channels recorded). 

(3) Spectral analyses are desired at a constant frequency resolution from zero fre- 
quency to the Nyquist frequency (half the sampling rate). 

(4) The filter warmup time cannot be tolerated. 

(5) The input data contain a large number of data samples and channels (large is on 
the order of 20 000 samples and 8 channels). 
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The filtering method is to be preferred when: 

(1) Correlation functions are not required for data reduction. 

(2) Cross -spectra are required on a large number of channels, 

(3) Spectral analyses are desired at a nonuniform frequency spacing (such as 
constant Q analysis or a high resolution analysis over a small section of the spectrum). 

(4) Auxiliary information (average frequency, equivalent bandwidth and equivalent 
duration in each frequency band) is required for stationary or randomness indicators. 

(5) Spectral answers are not desired at the Nyquist frequency and zero frequency. 

An additional difference in the spectral output between the two techniques is that 
the leakage effects due to side lobes of the estimation process in frequency can result in 
only positive power in the filtering approach, but in positive or negative power by the 
correlation technique. Therefore, leakage from adjacent frequencies can result in nega- 
tive coherences or high coherences (greater than 1) at low power levels in the correla- 
tion analysis, whereas the coherence function from the filtering analysis is always 
between 0 and 1. However, the filter approach may give a total integrated power that is 
greater than the total variance whereas the correlation method always insures that the 
integrated power is the total variance. 

Unless the user is reasonably certain about the frequency range of interest of his 
data, it is suggested that the correlation technique be used. This technique gives the 
user the maximum amount of information concerning the frequency content of his data. 
Thus, all sampled frequencies would be displayed, including any unexpected frequencies. 
This method uncovers possible errors in recording or digitizing the data. Also, the 
user may then specify certain frequency ranges for detailed analysis by using the fil- 
tering technique. 


Correlation Technique 


The correlation method also known as the Blackman-Tukey method is used to com- 
pute the spectral estimators. For a fixed data length T, the only spectral parameter not 
specified is the number of lags. 


The number of lags m is given by the following equation; 


m = 


1 

2 At Af 


(16) 


where Af is the desired equivalent resolution bandwidth (frequency resolution) for power 
spectral calculations. The term Af may have to be slightly altered so that m will be 
an integer. The maximum number of lags allowed by the TSA program is 999. See 
tables III and IV for the maximum number of channels to process at one time based on 
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the number of lags. In order to give good statistical reliability, it is desirable to keep 
the maximum number of la.gs m less than one-tenth the sample size. The sample size 
N is determined by T (time length of the data) and At by the relationship 

N = ~ (17) 

At 

The Blackman-Tukey method of computing the spectral estimators consists basically of 
computing the covariance function and performing a Fourier transform of this function. 
(See appendix A for the general computer procedure.) To eliminate the effect of taking, 
the transform of a truncated function, the covariance function (truncated in time at m At) 
is multiplied by a "lag window" as the transform is performed. This multiplication 
smooths the power spectral density, as explained in reference 3 (pp. 11 to 15). 

It should be noted that the TSA program does not distinguish between functions 
obtained from a time series with the mean subtracted and those obtained from a time 
series without the mean subtracted. Thus, one should exercise care when referring to 
other references on dynamic data about specific functions . For example, in the time 
domain functions , reference 2 denotes the function obtained from a time series with the 
mean subtracted as "autocovariance," the function obtained from a time series without 
the mean subtracted or with a mean of zero as "autocorrelation," and the normalized 
function with the mean subtracted as "Correlation coefficient." 

The equations for these functions and associated functions used by the TSA program 
are as follows: 

For autocovariance , 

N-k 

Rx(kAt)=j^ ^ (k= 0,1,2, . . ., m) (18) 

i=l 

where 

k the lag number 

m maximum lag number 

X channel 

k At displacement, r 

N number of data points 
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For autocorrelation, 

Rx(kAt) 

For cross -covariance, 


p (k At) 

Rx(0) 


/V _ A 1 9 


N-k 


Rxy(k At) - 2 


N - k Zo ^i^i+k 
i=l 


Rxy(“k At) = Ryjj(k At) — 


N-k 


yx' 


N - k ./ V 
1=1 


yi^i+k 


(k= 0,1,2, 


For cross-correlation. 


, , Rxv(±k At) 

p„,(±k at) = 'i/ - 2 

Ex(0)Ry(oy 


xy 


For power spectral density. 


ri 


N 


GxVl^'=2At 


m 

Rx(0) + 2 y XkRx(k At)cos(|^) 


k=l 


where 

r 

harmonic number 


frequency, f 

m 

X 

channel 

^k 

"lag window" (sm 


For Hamming smoothing, is defined as 

Tfk 


= 0,54 + 0.46 cos 


= 0.04 


m 


ml 1191 


* ? — / \ 


(k= 0, 1, 2, . . ., m) (20) 


m) (21) 


(k= 0, 1, 2, . . m) (22) 


(r = 0, 1, 2, . . ., m) (23) 


For Hanning smoothing, Xj^ is defined as 


Xi^= 0.50 + 0.50 cos 


Trk 

m 


(kss 1, 2, . . ., m-1) (24a) 
(k = m) (24b) 

(k = 1, 2, . , m-1) (25a) 
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Ar - 0 

For normalized power spectral density, 


g; 


/rf- 

Gx\-^ 


X \ m 


N 


Rx(0) 

For cross-spectral density 
/rf- 




31) =2 At 
m 


m 

D,,y(0) + 2 ^ At)COs(j^) 


k=l 


(k = m) 


(r = 0, 1, 2, . . ,, m) 
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where 


Dxy(k At) = ij 


D^y(k At) = i 


R^y(k At) + Ryx(k At) 
Rjjy(k At) - Ryx(k At) 


Gxy(f) = J[Cxy(fj] + [Qxy(f)]' 


Qvv(f) 

*xyW = 

For normalized cross-spectral density, 
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,N 


.N 


If C3x(f) or Gy(f) is negative, C^^(f) = Q^^(f) = 0. 


*xyW = ■“’xyW 


For coherence. 


r^Jl) = 


G^y(f) 


'’‘r' Gy(j) 

For amplitude transfer function 
Gxv(f) 


Txy(f) = 
Tyx(i) = ^ 


xy' 

Gx« 

Gx„(f) 


Gy(f) 


For inter power, 
m 
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(33) 

•> %) 

(34) 

M %) 

(35) 

M %) 

(36) 


(37) 


r=0 



(38) 


Filtering Technique 

The filtering technique involves the design of a band-pass filter for each frequency 
the user wishes to analyze. The data are passed through each filter and the output is 
squared and averaged to determine the spectral results. The only spectral parameters 
necessary for this process are the frequency specifications and the number of filter 
weights. 
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The program allows the following specifications for the frequency variable in the 

■f i Ho y'*! no' ‘-3 rin>'‘OCi ® 

^ JL^ ' X X i ^ ^ XX Vx JL £ e 

(1) Equispaced analysis; For equispaced analysis the frequencies are from fj to 
fmax 3-t equal intervals of Af. The filters will have bandwidths of Af and will be 

centered on fj, + Af, f j + 2 Af, . . and so forth, until (f| + m Af) ^ fmax- 


(2) Constant Q analysis: 
fj to fj^ax equal logarithmic intervals given as 
have center frequencies of f f j • 2^^^ 


For the constant Q analysis the frequencies are from 
k per octave. The filters will 


per octave . 
. vS/k 


and so forth, and 


bandwidths of f. 


2^/^ - 2“^/^y 2 , (2^/^ - l)/2 , (2^/^ - 2^/^) 


m/k 


and 


so forth, respectively, until f ^2 ' g ^max 


(3) Arbitrary analysis : For the arbitrary analysis, a table of frequencies 


f 


2 ’ 


fmax and the associated bandwidths Bj, B 2 , 


B 


max 


are read 


from the parameter deck. 


A restriction on the frequency variable is that f]^ cannot equal 0 and fmax can- 
not equal the Nyquist frequency. See tables V and VI for restrictions on the number of 
filters generated. Figures 6(a) to 6(c) illustrate each of these band-pass filter options. 
The plots shown are used for descriptive purposes only and are not response functions of 
the actual analysis filters. 

The analysis filters are cascades of simple band-pass filters. The outputs are 
computed by the use of recursive equations on the past outputs; therefore, they have an 
initial settling time when the data are "turned on.” The program allows the number of 
data samples equal to the reciprocal of the normalized filter bandwidth for a warmup 
period before the filter outputs are used in spectral computations. 

For cross-spectral computations, both "in-phase" (real) and "quadrature" (imagi- 
nary) filters are required. (See fig, 7.) These filters are produced by including filters 
in the cascade that shift the phase +45° and -45° at the center frequency. The center fre- 
quency of the filter is shifted to the low-frequency side for the in-phase and to the high- 
frequency side for quadrature phase. By setting the proper input parameters, the pro- 
gram will print and/or plot the desired filter response functions. 

The program allows the user a choice of the number of filter weights used to con- 
struct the required band-pass filters. The user can choose between 5, 7, or 9 weights 
which corresponds to a cascade of 2, 3, or 4 filters. 

Naturally, a sharper and more precise filter can be constructed when using more 
filters in the cascade, but this procedure requires more computer time and storage. 
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Therefore, a trade-off must be made in terms of filter characteristics and computer 
time and storage. See tables ¥ and ¥I for restrictions on the number of filter weights. 

Each of the filters has the following recursive weights: 


bj(fo) = 2 exp(-7ry)cos<^27ro>Q 
^2(^0) = -exp(-27Ty) 


Z_ 

2 CO, 


(39a) 

(39b) 


where 




0.4f' 


y 


B 




1/2 


0.8f’ 


(Uo = < 1.2f 




(For first filter in cascade) 
(All remaining filters in cascade) 

(For first filter in cascade (in-phase)) 
(For first filter in cascade (quadrature)) 
(All remaining filters in cascade) 


f 

Vr 


normalized filter bandwidth associated with Iq 
number of filters in cascade 

R/ 


For both the in-phase filter recursive weights Bj (fo), and quadrature fi 


filter recursive 


\th r-t 


weights B^(fQ), the ith cascaded filter is built up from the (i - 1) filter by the following 
relationship: 


BfK) = - bi(y B|T)(fo) - b2(fo) 

with the initial conditions : 

sj/g) = 0 


(i = 1, 2, , . 2i) (40) 


(i 0, any i) 
(Any i) 


The lone nonrecursive weight Ao(fo) for the cascaded filter for both the in-phase and 
quadrature filters is computed by 
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where 


ko(io) 


¥-^l 

1 - y Bi(fQ)exp[-ij2iTfQ] 

j=l 


(41) 


V number of filter weights 

i = \pr 


The application of the band-pass filters is then 



V-1 

Ao('o)]=‘n+ ^ [Br(fo)|k.i(fo) 

i=l^ 

(n = 1 , 2 , . 

. N) 

(42) 


V-1 

Ao('°) \ Pi('oK-i(y 

i=l 

(n = 1, 2, . 

. ., N) 

(43) 


In order to minimize computer processing time and reduce disk storage, the program 
reduces the sampling rate of the data as much as possible for each filter output. These 
filter decimation factors are checked to insure that the new equivalent sampling rate is 
not a multiple of the filter's center frequency. Also, there must be, after reduction, a 
minimum of 50 points to analyze and 5 points per center frequency cycle. The filter 
decimation factors D(fo) are determined in the following manner: 

Initially, 



where 


dp precondition decimation factor 

fg normalized filter bandwidth associated with fg 

I(x) integer part of x 

Define 

M(fo) = l(2f;D(fo) dp) 

where 


normalized filter center frequency 


( 45 ) 
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If M(fo) is even, 


Xo = foD(fo) dp - i M(fo) 

If M(fo) is odd, 

Xo = |(M(fo) + l) - f;D(fo) 


dr 


(47) 


If the following conditions are not satisfied, decrease D(fo) by 1, recompute M(fo) and 
Xq, and try again: 

5 ^ . (48) 


Xo + ffB^O.5 


X, 


5 r' 


O 4 "B 

N 


f^ a 0 


D(fo) 


50 


(49) 


(50) 


where N-^ = II 


1 b dpj 


Number of filter warmup samples, 


^ 5 


fo D(fo) 

(fg/f) modulo D(fQ) # 0 


(51) 

(52) 


The conditions are checked in the order that they appear. If D(fo) = 1 and either 
condition (48) or (49) fails, a diagnostic message is printed and the job is aborted. If 
D(fo) = 1 and either condition (50) , (51), or (52) fails, a diagnostic message is printed and 
the program continues to execute. 

When the spectra are calculated by the band-pass filtering technique, there are no 
correlation computations performed, but additional information is obtained in the form of 
averages on the filter outputs (equivalent bandwidth and equivalent duration). 

The equivalent bandwidth is a measure of the bandwidth of the filter as determined 
by the filter output. The bandwidth of each filter is computed as the root mean square of 
the instantaneous frequency around the average frequency. Small equivalent bandwidths 
(on the order of one-tenth of the filter bandwidth or less) may be the result of poor pre- 
whitening, or the presence of a sinusoid. 
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The equivalent duration is a measure of the fraction of the time that the output of a 

filter is at high power. It is defined as the ratio of average power to peak power output 
of a filter. Thus , a sine wave would have an equivalent duration of 0.5, random noise 
would be expected to be in the neighborhood of 0,1, and transients or nonstationary da.ta 
would be much smaller. 


The equations for the spectral computations (noted only when different from the 
correlation technique) are as follows: 


Power spectral density: 
M _ 


^x(fo) “ 


I K('c)) - 
n=l^ 

2M|/3(fo) 


where 


(53) 


M number of filter outputs after decimation 


^(fo) 


noise bandwidth of the filter pair (computation method shown in appendix B) 


Cross-spectral density: 

M 

I \K(^o) -H iSy iSfo) 

Cxy(^o) = 


M^(fo) 


(54) 


M 


I^'o) - ^ (fo) RSto) 


^xy(^o) 


n=l 




(55) 


In the equations for average frequency and equivalent bandwidth below, A$^(fo) is 
defined by 




i*hfo) - 


tan" 


tan 


In both cases, is shifted to lie between 






tan" 


-■nPo)^ 


K('o), 


rX 


tan 


■l/ "^n-l(^oc 
vRn-l(to); 


27rAr 


2t7Ao 


^M(fQ) even^ 


M(fo) odd 


) ( 56 ) 


J 


TT and 'IT. 



Average frequency: 


A* A A 

(57) 

where 


e(fo) = 2tt AtMD(fo) dp 


M 






Equivalent bandwidth: 


' |4,2(At) W(y d2(M - 1) 

(58) 

where 


M , 

1 (^^n(«) 

A#2(f ) = 2 zL_ — 

^ ' M 


Equivalent duration: 


„ ^ (5(to) Gx(fo) 

Dx(fo) = p. ^ ^ 

_^n(^o)]max 

(59) 

where 


<((o)= KA)^(*o))' 


Filter response function: 


r 2 2 

(60) 

where the range of f is given by the parameter cards and 

oq 

Cj kJ 





A qPo) 


a R (in-phase) or I (quadrature) 

i= fl 


$j^(f) = tan“l 



imag 


real 


A#f^(f’) = #f^(f ) - #J^(f) 


MODIFICATION OF SPECTRAL OUTPUTS 


Several options are available for modification of the spectra outputs to reflect 
more accurately the spectra of the data of interest. These options may be chosen inde- 
pendently of each other and of the method chosen for obtaining the spectra. In each of 
the options where the power spectra and cross -spectra are modified, the functions which 
are determined by these spectra are computed by use of the modified spectra. The 
options are listed and the spectra modification equations are given in this section. 


Amplitude and Phase Calibration Tables 


A calibration of the instrumentation system amplitude response H(f) and/or phase 
response j3(f) may be entered as a table of the variation of response with frequency for 
each channel. If the frequencies in the calibration tables do not match the frequency 
needed for the spectra, a linear interpolation is performed to obtain the proper response 
for the spectral frequency. This information is used to correct all spectral outputs for 
each channel in the following manner; 


Gx(f) = 




(61) 


Cxy(f) 


Cxy(f) 

Hx(f) Hy(f) 


(62) 
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( 63 ) 


o (f) 

'xy ■ ■ Hx(f) Hy(f) 


( 64 ) 


Time Delay Calibration 

A sinusoidal calibration signal may be applied simultaneously to all channels of the 
input data at time of data recording. If this signal is present, it may be used to compute 
a time lag r^y between each pair of channels x and y for which cross-spectra or cross- 
correlations are desired. The first 10 (or less if there are not 10) zero crossings are 
found by interpolation; they are then averaged and the difference between each set of 
two channels is taken as the time lag between those two channels. This time lag is then 
used to correct the phase of cross -spectral functions as follows: 

^xyW = ^xy(^) - 360fTxy (65) 

Generally, the time delay calibration is not necessary unless the planned digitizing 
sampling frequency is very high and hence makes the possible small time difference 
between channels significant. 


Shaping Filter Postdarkening 

If the data were prewhitened with a shaping filter, the effect of this filter may then 
be removed from the spectra by using the transfer function of the filter T(f). If the 
data were treated with a shaping filter to remove undesirable frequencies, the spectra 
should not be modified by the transfer function of the filter since this modification would 
simply insert the effect of these frequencies on the spectra, as if they were not removed. 
The filter transfer function is determined by 


2 _ PnW ^ OnW 

P,?(f)*Q2(f) 


where 


N 

y 

/- 7 

i=o 


Pj^(f) = y Aj cos(27rjf) 


Ij^(f) = 2 ^ Aj sin( 27 rjf' 

3=0 


( 66 ) 


(67) 


( 68 ) 
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(69) 


D 

Pj^(f) = 1 - Bjii cos(27mif*) 
m=l 

D 

Qj^(f) = Bjjj sin(2Timf') 
m=l 

If desired, the spectra are postdarkened as follows: 


Gx(f) = 


G^(f) 

T^(f) 


Cxy(f) 


T^(f) 


Qxy(f) = 

^ T2f) 


PROBABILITY DENSITY ESTIMATION METHODS 


(70) 


(71) 


(72) 


(73) 


Two techniques are available for the estimation of probability densities: (1) the 
standard histogram technique of dividing the amplitude of the data into intervals (bins) 
and counting the number of data points that fall in each bin, and (2) the smooth density 
function estimate obtained by an orthogonal series expansion. 

The primary defects of the histogram method are: 

(1) The division of the data into arbitrary bins may distort the distribution if there 
is curvature within a bin. 

(2) A large number of bins are required to define the shape of a distribution and, 
as a result, a very jagged histogram is obtained. 

The orthogonal function method does not have these defects but 

(1) The distribution must be smooth to converge with 20 coefficients (orthogonal 
functions) . 

(2) The general shape of the distribution must be known. 

(3) The approximation and "leakage" effects may cause some probability values to 
be negative. 

(4) The chi-square normality test cannot be performed. 
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Unless the user is reasonably certain about the shape of the probability density 
function of his data, it is suggested that the histogram technique be used. This technique 
gives the user the shape of the probability density function, and the orthogonal function 
technique could then be used for obtaining a smooth probability distribution, if desired. 

Histogram Technique 

A histogram is computed for a channel of data by dividing the amplitude of the data 
into a specified number of bins. The number of data samples that lie in each bin are 
counted, and a plot is produced with data amplitude on the X-axis and counts on the Y-axis. 
If the parameter to subtract the mean was selected, then the mean was subtracted before 
the histogram was computed. Hence, to convert the histogram amplitude to actual data 
amplitude, one must add the mean for that channel. Provisions are made to combine 
adjacent bins in order to reach a specified number of counts per bin. This procedure 
tends to smooth the histogram at the low-count regions. Figure 8 gives an illustration 
of some digital input data and its corresponding histogram. 

There are four options available for specifying the histogram parameters. For all 
four options , the number of bins is an input quantity. For option 1, the range of the histo- 
gram is determined by specifying the maximum and minimum histogram amplitude values; 
for option 2, the program scans the data for the maximum and minimum values and uses 
these values to determine the range of the histogram; for option 3, the maximum and 
minimum values are specified and bins are combined to reach an input-specified minimum 
number of counts in a bin; for option 4, the program scans the data for the maximum and 
minimum values and bins are combined as in option 3. 

There is a restriction on the total number of bins for the input channels. If the 
correlation technique is used for spectral analysis, the total allowable number of bins 
is 1000. If the filtering technique is used for spectral analysis, the total number of bins 
allowed is reduced to 800. 


Orthogonal Function Technique 

By expanding the probability density into a finite orthonormal series and estimating 
the orthonormal coefficients from the data, a smooth estimate of the probability density 
function can be produced. The following two expansion choices are offered; Hermite 
functions which are appropriate for two-tailed distributions such as the Gaussian, chi- 
square for large degrees of freedom, and so forth; Laguerre functions which are appro- 
priate for single-tailed distributions such as the exponential, gamma, or chi-square dis- 
tribution for small degrees of freedom. 
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The equations for computing the orthonormal Hermite functions are as follows; 


exx)(-z?l^ 


7T 


1/4 


Fi(Zi) = ^Zj^FQ(Zi) 


( 75 ) 


where 


Fj(^i) = f - (V-J 0 = 2,3,..., 


= 


X. -X 


19) 


(76) 


The plots of these functions are shown in figure 9. The equations for computing the 
orthonormal Laguerre functions are as follows: 

/-zh 

Fo(Zi) = expj^— 

Fi(Zi) = (1 - Zi)FQ(Zi) 

(2j - 1 - Zi)Fj_i(Zi) - (j - l)Fj_2(Zi) 

Fj(Zi) = 1 J L£L_i (j = 2, 3, . . ., 19) 


(77) 

(78) 

(79) 


where 


Zi = 


_ 2(xj - b) 
X - b 


r minimum Xj^ if tail goes toward maximum 


(jmaximum Xj if tail goes toward minimum 


The plots of these functions are shown in figure 10. The coefficients for either expansion 
series are computed by 


N 


N 


F](Zi) 


(80) 


i=l 
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whei-e 


N the number of data points 

normalized data point 

The normalized probability density function P(z) over the amplitude range specified in 
the input is then computed by a series expansion by using the appropriate series and the 
coefficients determined. The equation for this function is 

P 

P(Zi) = ^ ajFj(Zi) (81) 

3=1 

where 

P number of expansion terms (input parameter) 

Zj normalized Xj 

z determined by the range and increment specified on input cards 

The probability density function P(x) is then obtained from the normalized probability 
density function in the following manner: 

For Hermite functions , 

P(x) = ^ (82) 

For Laguerre functions, 

P(x) = (83) 

|X - b| 

Note that P(x) and P(z) are computed so that the area under each curve equals 1. 

ADDITIONAL OUTPUTS 

Many additional or subsidiary outputs are provided by the program for checking the 
form of the probability function, the digitizing process, off-scale values, randomness, and 
stationarity. 
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Averse on the Entire Data 


The mean, standard deviation, third (skewness), and fourth (kurtosis) normalized 
central moments of the data are computed by selecting the proper input parameter values. 
These values are always computed on the raw input data before preconditioning, if 
required, was performed. The equations for these values are as follows: 

Mean: 

N 

X = (84) 

i=l 

Standard deviation: 


a = 



Skewness: 


/3 = 



Kurtosis: 


K = 



(85) 


( 86 ) 


(87) 


An interpretation that can be applied to the third and fourth normalized central moments 
is found in the discussion on "Distribution Tests." 


Cross-Correlation Time Lag 

If a time delay calibration (discussed earlier) is present, the true maximum time 
lag for cross-correlations may be computed. This value is based upon the time lag- 
found between the two channels and the normal maximum time lag. Normally, the maxi- 
mum time lag for cross-correlations is found by simply multiplying the maximum number 
of lags by the time increment. Therefore, the true maximum time lag between two chan- 
nels is determined as follows: 

= m At - T (88) 

xy xy ^ ' 
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where 


m maximum number of lags 


Tj^y time lag between channels x and y 


Distribution Tests 

Several of the program's outputs may be used to test the data for the general shape 
of its distribution. For a two-tailed distribution, the shape is generally compared with 
the Gaussian distribution. For a single-tailed distribution, the shape is generally com- 
pared with the exponential distribution. 

The third and fourth normalized central moments can give a rough idea of the dis- 
tribution. To give the user some insight into the meaning of these moments, a Gaussian 
distribution has a skewness value of 0 and a kurtosis value of 3. A negative skewness 
value indicates the distribution has an elongated left tail and a positive skewness value 
indicates an elongated right tail. A distribution is called mesokurtic, platykurtic, or 
leptokurtic as its kurtosis equals 3, less than 3, or greater than 3, respectively. (See 
ref. 4.) A platykurtic distribution indicates the data points are less heavily concentrated 
about the mean and a leptokurtic distribution indicates the data points are more heavily 
concentrated, as compared with the Gaussian distribution. 

The chi-square goodness of fit to a Gaussian distribution may be computed when 
the histogram option is chosen for the probability density function. The method of com- 
puting the chi-square value y to test for normality is as follows: 


where 


z — 

i=l 


(89) 


probability of being in ith bin if Gaussian 
Nj^ number of points in ith bin 

N total number of points 


K degrees of freedom, number of bins - 2 


The region of normality acceptance is 



(90) 

QQ 



2 

where the ¥alue of X|^ o, is available from table ¥II and a is the proba,bility of a, 

type 1 error. As an example; suppose y = 65; K = 50; and the data are to be tested for 

normality at the a = 0.05 level of significance; then the data, would be accepted to be 
2 

norma,l since Xi^q q5 ” 67.5, (See ref. 5.) 

When using the orthogonal function option for the probability density function, one 
can use the orthogonal coefficients to estimate the normality of the data. If the data are 
Gaussian, the zero-order coefficient in the Hermite expansion is three or four orders of 
magnitude greater than the remaining coefficients. If the data are exponentially dis- 
tributed, the zero-order coefficient in the Laguerre expansion is three or four orders of 
magnitude greater than the remaining coefficients. 

Accuracy Test 

If the spectra were computed by the correlation technique, the program has an option 
for computing an accuracy measurement for the spectral estimators. For a specified 
confidence percentage, the program computes a percentage band for the spectral estima- 
tion. For example, a confidence percentage of 98 and a percentage band of 4.0 means that 
the user can be 98 percent certain that the spectral answers are within 4.0 percent of the 
true value. The equation for the percentage confidence band is as follows; 



where 

N number of data points 

m maximum lag number 

decibel spread for specified confidence percentage (see table VIII) 

Off-Scale Test 

The user can input by means of parameter cards an off-scale value and a maximum 
number of off-scale points for a channel. The program then checks each data point to 
see whether its absolute value exceeds the off-scale value . If it does, the data point , 
time, and replacement value is printed. When the maximum number of off-scale points 
for one channel is exceeded, the program aborts that case. The replacement value is the 
mean of all the preceding points. If the off-scale point is the first point, the replacement 
value equals zero. 
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The maximum and minimum value for each channel are determined after all off- 
scale values have been replaced and, if requested, the bias removed. These values are 
automatically printed to give the user the range of the data analyzed. 

CONCLUDING REMARKS 

This paper was prepared as a reference guide to users of the Langley time series 
analysis computer program (designated R1299), which is a general purpose program for 
the analysis of random, stationary time series. The bases for choosing available options, 
descriptions of the options, possible output quantities, equations of the various computa- 
tions, and tables and figures which might aid in the analysis of data were included. 

Langley Research Center, 

National Aeronautics and Space Administration, 

Hampton, Va., December 1, 1970. 



APPENDIX A 


GENERAL COMPUTER PROCEDURE 

The program is set up so that any of the 20 channels may be processed according 
to parameter card options. Cross-correlations or cross-spectra may be computed 
between all possible combinations of any two channels. 

The following procedure is used: 

(a) All parameters and calibration data are read and the time lag, if any, is 
computed. 

(b) The data file is read and the mean, if requested, is computed for each channel. 

The data are placed on a blocked disk file. 

(c) When phase II (filtering method) is used, the filter weights are computed and the 
filter response, if requested, is computed, printed, and placed on a binary output tape. 

(d) A data point is read from the data disk file. 

(e) Probability density function estimators, if requested, are updated for all channels. 

(f) The statistics (standard deviation, third central moment, and so forth) summa- 
tions for all channels are updated. 

(g) The data point is conditioned according to parameter card options. 

(h) This step is the first one in which all the channels may not be processed 
together. The number of channels processed in one data pass is a function of the num- 
ber of lags in the correlation method and of the number of weights and filters in the 
filtering method. 

For correlation method: 

(1) The autocovariance and cross - covariance sums for channels to process 

during this data pass are updated. 

(2) If all the data have not been read, return to step (d). 

For filtering method: 

(1) The channels to process during this data pass are filtered and placed 

on a scratch disk file. 

(2) If all the data have not been read, return to step (d), 

(i) After all the data have been read, the probability density function estimations, 

if requested, are completed, printed, and placed on the binary output tape. The statistics 
are completed and printed. 
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APPENDIX A - Continued 


(j) For correlation method: 

(1) The autocovariance anc 

printed, and placed on the binary output tape. 


ins a] 


)rfipl0t6(3 ^ 


(2) The autocorrelation, cross-correlation, spectra, etc., are computed, 

printed, and placed on the binary output tape. 

(3) If there are more channels to process (that is, another data pass is 
necessary) , rewind the data disk file and return to step (d) and eliminate steps (e) , 
(f) , and (i) . 


For filtering method: 

(1) If there are more channels to process (that is, another data pass is 
necessary), rewind the data disk file and return to step (d) and eliminate steps (e), 
(f), and (i). 

(k) Filtering method only; 

(1) The spectral summations are computed from the filtered data on the proper 
scratch disk file. 


(2) The final spectra, equivalent bandwidth, and so forth answers are com- 
puted, printed, and placed on the binary output file. 

(3) Return to step k(l) if more than one data pass was required to filter all 
the channels. 


(1) After all the channels have been processed, the plot subroutines are entered 
and the plot instructions are generated by use of the binary output file. If another data 
case is present, the program will return to step (a); if not, it will terminate. 

For the correlation method, the number of channels processed in one data pass in 
steps (h) and (j) depends upon the maximum number of lags desired. There are 2000 
storage locations allotted for single -channel processing and 2000 storage locations 
allotted for cross-channel processing. Tables III and IV show the maximum number of 
single channels and cross channels that can be performed in one data pass for a given 
number of lags. 

For the filtering method, the number of channels processed in one data pass in 
steps (h), (j), and (k) depends upon the number of filter weights and the number of filters 
that must be generated for the frequency analysis. There are 2000 storage locations 
allotted for the filter weights and 2000 storage locations allotted for storing the previous 
filter outputs from each channel. Tables V and VI show the maximum number of single 
channels or cross channels that can be performed in one data pass for a given number of 
filter weights and filters that must be generated. When possible, an even number of 
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APPENDIX B 


NOISE BANDWIDTH COMPUTATIONAL METHOD FOR FILTER PAIR 


The equation for computing the noise bandwidth for the analysis pair of filters used 

in the filtering technique is as follows; 

■ 1/2 


/3(fo) 


At dp oq 
1 


At di 


P 


H(fQ,G')dQ! 

Q!o(fo)-Aa!(fo) 


0 


i 


1/2 

Q!o(fo)+Aa!(fo) 


H(fo,ff)do! + 
H(fo,«)dQ! 


‘“o(y+Ao!(fo) 

Q'o(fo)-AQ'(fo) 


H(fo,Q!)dQ' 


where 



«o(fo) = fo At dp 



k(Bi) 

(fo = «l) 

A«(fo) = ^ 

1 

i^o(^o) " '^o(^o) (^0 “ ^2’ ^3» • • 


Bl 

filter bandwidth associated with fj 


fl 

center frequency of first filter 


f2 

center frequency of second filter 


fo 

center frequency of filter preceding Iq 



P(fo,0!) +Q(fo,0!) 
H(fo,«) = — — 


P(fo,Q') 


A IR./ r V 

o{^o) 


i y~i 


/V=l 


B^(fo) cos(27rai) | + ( ^ B"(fo) sin(2ira'i) 


i=l 


, i=l 
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APPENDIX B - Concluded 



Simpson's rule for integration with 90 intervals is used in each of the 3 regions to evalu- 
ate the integrals and obtain the final result. 
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TABLE L- DYNAMC DATA RECORD FORMAT FOR TSA PROGRAM 


Word number 

Mode^ 

Description 

1 

Floating 

Number of channels in the record 

2 

Floating 

Serial number 

3 

Fixed (Hollerith)^ 

Primary engineering identification 

4 

Fixed (Hollerith)^ 

Primary engineering identification 

5 

Fixed (Hollerith)*^ 

Additional engineering identification 

6 

Fixed (Hollerith)^ 

Additional engineering identification 

7 

Fixed (Hollerith)*^ 

Additional engineering identification 

8 

Fixed (Hollerith)^ 

Additional engineering identification 

9 

Floating 

Record count 

10 

Floating 

Elapsed time in seconds 

11 

Floating 

Channel 1 data sample 

12 

Floating 

Channel 2 data sample 

N + 10 

Floating 

Channel N data sample where N is 
the number in word 1 


^Data should be written on tape as a binary file using a floating- 


point array. 

^Each word contains 10 characters. 


TABLE n.- BLOCKED FORMAT FOR TSA PROGRAM 


Number of 
channels, 
N 

Words per 
record 

Records per 
block^^ 

Words per 
block^ 

o 

VII 

VII 

1—1 

20 

25 

500 

10 < N i 20 

30 

17 

510 


block is a FORTRAN logical record. 





TABLE III.- MAXIMUM NUMBER OF CROSSES PROCESSED 
IN ONE PASS FOR CORRELATION METHOD 


Number of 
lags desired 

Maximum number of crosses 
processed in one data pass 

500 to 999 

1 

334 to 499 

2 

250 to 333 

3 

200 to 249 

4 

167 to 199 

5 

143 to 166 

6 

125 to 142 

7 

112 to 124 

8 

100 to 112 

9 

91 to 99 

10 




TABLE IV.- MAXIMUM NUMBER OF SINGLE CHANNELS PROCESSED 
IN ONE PASS FOR CORRELATION METHOD 


Number of 
lags desired 

Maximum number of single channels 
processed in one data pass 

667 to 999 

2 

500 to 666 

3 

400 to 499 

4 

334 to 399 

5 

286 to 333 

6 

250 to 285 

7 

223 to 249 

8 

200 to 222 

9 

182 to 199 

10 

167 to 181 

11 

154 to 166 

12 

143 to 153 

13 

134 to 142 

14 

125 to 133 

15 

118 to 124 

16 

112 to 117 

17 

106 to 111 

18 

100 to 105 

19 

0 to 99 

20 









t TABLE VI.- MAXIMUM NUMBER OF SINGLE CHANNELS PROCESSED 
IN ONE PASS FOR FILTERING METHOD 


Number of 
filtering weights 

Number of filters 
generated 

Maximum number of 
single channels processed 
in one data pass 

5 

201 to 400 

1 

5 

101 to 200 

2 

5 

67 to 100 

4 

5 

51 to 66 

6 

5 

41 to 50 

8 

7 

143 to 285 

1 

7 

72 to 142 

2 

7 

48 to 71 

4 

7 

36 to 47 

6 

7 

29 to 35 

8 

9 

112 to 222 

1 

9 

56 to 111 

2 

9 

38 to 55 

4 

9 

28 to 37 

6 

9 

23 to 27 

8 
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GOODNE SS “O F* =" ^'1 T TE S T 


K 

a = 0.20 

a = 0.15 

a = 0.10 

a = 0.05 

a = 0.02 

a = 0.01 

a = 0.005 

a = 0,001 

a = 0.0005 

1 

1.64 

2.03 

2.71 

3.84 

5.41 

6.64 

7.88 

11.2 

12.6 

2 

3.22 

3.76 

4.60 

5.99 

7.82 

9.21 

10.60 

14.1 

15.7 

3 

4.64 

5.29 

6,25 

7.82 

9.84 

11.3 

12.80 

16.6 

18.1 

4 

5.99 

6.72 

7.78 

9.49 

11.7 

13.3 

14.9 

18.7 

20.4 

5 

7.29 

8.09 

9.24 

11.1 

13.4 

15.1 

16.8 

20.8 

22.4 

6 

8.56 

9.42 

10.6 

12.6 

15.0 

16.8 

18.6 

22.7 

24.4 

7 

9.80 

10.7 

12.0 

14.0 

16.6 

18.5 

20.3 

24.5 

26.3 

8 

11.0 

12.0 

13.4 

15.5 

18.2 

20.1 

22.0 

26.3 

28.1 

9 

12.2 

13.3 

14.7 

16,9 

19.7 

21.7 

23.6 

28.1 

29.9 

10 

13.4 

14.5 

16.0 

18.3 

21.2 

23.2 

25.2 

29.8 

31.7 

11 

14.6 

15.8 

17.3 

19.7 

22.6 

24.8 

26.8 

31.4 

33.4 

12 

15.8 

16.9 

18.5 

21.0 

24.1 

26.3 

28.3 

33.1 

35.0 

13 

17.0 

18.2 

19.8 

22.4 

25.5 

27.7 

29.8 

34.7 

36.7 

14 

18.2 

19.4 

21.0 

23.7 

26.9 

29.1 

31.3 

36.3 

38.3 

15 

19.3 

20.6 

22.3 

25.0 

28.3 

30.6 

32.8 

37.8 

39.9 

16 

20.5 

21.8 

23.5 

26.3 

29.6 

32,0 

34.3 

39.4 

41.5 

17 

21.6 

23.0 

24.8 

27.6 

31.0 

33.4 

35.7 

40.9 

43.1 

18 

22.8 

24.1 

26.0 

28.9 

32.3 

34.8 

37.2 

42.5 

44,6 

19 

23.9 

25.3 

27.2 

30.1 

33.7 

36.2 

38,6 

44.0 

46.2 

20 

25.0 

26.5 

28.4 

31.4 

35.0 

37.6 

40.0 

45.4 

47.7 

21 

26.2 

27.7 

29.6 

32.7 

36.3 

38.9 

41.4 

46.9 

49.2 

22 

27.3 

28.8 

30.8 

33.9 

37.7 

40.3 

42.8 

48.4 

50.7 

23 

28.4 

30.0 

32.0 

35.2 

39.0 

41.7 

44.2 

49.9 

52.2 

24 

29.6 

31.1 

33.2 

36,4 

40.3 

43.0 

45.6 

51.3 

53.7 

25 

30.7 

32.3 

34.4 

37.6 

41.6 

44.3 

46.9 

52.7 

55.1 

26 

31.8 

33.4 

35.6 

38.9 

42.9 

45.6 

48.3 

54.2 

56.6 

27 

32.9 

34.6 

36.7 

40.1 

44.1 

47.0 

49.6 

55.6 

58.0 

28 

34.0 

35.7 

37.9 

41.3 

45.4 

48.3 

51.0 

57.0 

59.4 

29 

35.1 

36.8 

39.1 

42.6 

46.7 

49.6 

52.3 

58.4 

60.9 

30 

36.3 

38.0 

40.3 

43.8 

48.0 

50.9 

53.7 

59.8 

62.3 

31 

37.4 

39.1 

41.4 

45.0 

49.2 

52.2 

55.1 

61.2 

63.7 

32 

38,5 

40.3 

42.6 

46.2 

50.5 

53.5 

56.4 

62.6 

65.1 

33 

39.6 

41.4 

43.7 

47,4 

51.8 

54.8 

57.7 

64.0 

66.5 

34 

40.7 

42.5 

44.9 

48.6 

53.0 

56.1 

59.0 

65.4 

67.9 

35 

41.8 

43.6 

46.1 

49.8 

54.3 

57.4 

60.3 

66.7 

69,3 

36 

42.9 

44.8 

47.2 

51.0 

55.5 

58.6 

61.6 

68.1 

70.7 

37 

44.0 

45.9 

48 »4 

52.2 

56,7 

59.9 

62.9 

69.4 

72.1 

38 

45.1 

47.0 

49.5 

53.4 

58.0 

61.2 

64.2 

70.8 

73.5 

39 

46.2 

48.1 

50.7 

54.6 

59.2 

62.5 

65.5 

72.2 

74.9 

40 

47.3 

49.2 

51.8 

55.8 

60.4 

63.7 

66.8 

73,5 

76.2 









TAtsL 


K 

a = 0.20 

O' = 0.15 

41 

48.4 

50.4 

42 

49.5 

51,5 

43 

50.5 

52.6 

44 

51.6 

53.7 

45 

52.7 

54.8 

46 

53.8 

55.9 

47 

54.9 

57.0 

48 

56.0 

58.1 

49 

57.1 

59.2 

50 

58.2 

60.3 

51 

59.2 

61.4 

52 

60.3 

62.5 

53 

61.4 

63.7 

54 

62.5 

64.8 

55 

63,6 

65.9 

56 

64.7 

66.9 

57 

65.7 

68.0 

58 

66.8 

69.1 

59 

67.9 

70.2 

60 

69.0 

71.3 

61 

70.0 

72.4 

62 

71.1 

73.5 

63 

72.2 

74.6 

64 

73.3 

75.7 

65 

74.3 

76.8 

66 

75.4 

77.9 

67 

76.5 

79.0 

68 

77.6 

80.1 

69 

78.6 

81.2 

70 

79.7 

82.3 

71 

80.8 

83.3 

72 

81.9 

84.4 

73 

82.9 

85.5 

74 

84.0 

86.6 

75 

85.1 

87.7 

76 

86.1 

88.8 

77 

87.2 

89.9 

78 

88.3 

90.9 

79 

89.3 

92.0 

80 

90.4 

93^1 


E viL- CHI-tQUARE ( 


a = 0.10 

a = 0.05 

52.9 

56.9 

54.1 

58.1 

55.2 

59.3 

56.4 

60.5 

57.5 

61.7 

58.6 

62.8 

59.8 

64.0 

60.9 

65.2 

62.0 

66.3 

63.2 

67.5 

64.3 

68.7 

65.4 

69.8 

66.5 

71.0 

67.7 

72.2 

68.8 

73.3 

69.9 

74.5 

71.0 

75.6 

72.2 

76.8 

73.3 

78.0 

74.4 

79.1 

75.5 

80.2 

76.6 

81.4 

77.7 

82.5 

78.9 

83.7 

80.0 

84.8 

81.1 

86.0 

82.2 

87.1 

83.3 

88.3 

84.4 

89.4 

85.5 

90.6 

86.6 

91.7 

87.7 

92.8 

88.8 

93.9 

90.0 

95.1 

91.1 

96.2 

92.2 

97.4 

93.3 

98.5 

94.4 

99.6 

95.5 

100.8 

96.6 

101.9 


i;S'r - Continued 


a = 0.02 

a = 0.01 

a = 0.005 

a = 0.001 

a = 0.0005 

61.7 

65.0 

68.1 

74.8 

77.6 

62.9 

66.2 

69.4 

76.2 

78.9 

64.1 

67.5 

70.7 

77.5 

80.3 

65.3 

68.7 

71.9 

78.8 

81.6 

66.6 

70.0 

73.2 

80.2 

83.0 

67.8 

71.2 

74.5 

81.5 

84.3 

69.0 

72.5 

75.7 

82.8 

85.7 

70.2 

73.7 

77.0 

84.1 

87.0 

71.4 

74.9 

78.3 

85.4 

88.3 

72.6 

76.2 

79.5 

86.7 

89.7 

73.8 

77.4 

80.8 

88.1 

91.0 

75.0 

78.6 

82.0 

89.4 

92.3 

76.2 

79.9 

83.3 

90.6 

93.6 

77.4 

81.1 

84.5 

92.0 

94.9 

78.6 

82.3 

85.8 

93.3 

96.3 

79.8 

83.5 

87.0 

94.5 

97.6 

81.0 

84.8 

88.3 

95.8 

98.9 

82.2 

86.0 

89.5 

97.1 

100.2 

83.4 

87.2 

90.8 

98.4 

101.5 

84.6 

88.4 

92.0 

99.7 

102.8 

85.8 

89.6 

93.2 

101,0 

104.1 

87.0 

90.8 

94.5 

102.2 

105.4 

88.1 

92.0 

95.7 

103.5 

106.7 

89.3 

93.2 

96.9 

104.8 

108.0 

90.5 

94.4 

98.1 

106.1 

109.3 

92.7 

95.6 

99.4 

107.3 

110.5 

92.9 

96.9 

100.6 

108.6 

111.8 

94.0 

98.1 

101.8 

109.9 

113.2 

95.2 

99.3 

103.0 

111.2 

114.4 

96.4 

100.4 

104.3 

112.4 

115.7 

97.6 

101.6 

105.5 

113.7 

116.9 

98.7 

102.8 

106.7 

114.9 

118.2 

99.9 

104.0 

107.9 

116.2 

119.5 

101.1 

105.2 

109.1 

117.4 

120.8 

102.3 

106.4 

110.3 

118.7 

122.0 

103.4 

107.6 

111.5 

119.9 

123.3 

104.6 

108.8 

112.7 

121.2 

124.6 

105.8 

110.0 

113.9 

122.4 

125.8 

106.9 

111.2 

115.2 

123.7 

127.1 

108.1 

112.4 

116.4 

124.9 

128.4 











i/vrsLiH. vii.“ 


A T-iT71 y~«^/^T-VXTT?lC«0 Tt1 TTrr rT» 

“i^uAnjii Lj\-/<^i-/r<i!iOO“V^r “rl X xjiroj 


A^uiiCiuueu 


K 

a = 0.20 

a = 0.15 

0 ! = 0.10 

a = 0.05 

a = 0.02 

a = 0.01 

ffi = 0.005 

a = 0.001 

a = 0.0005 

81 

91.5 

94,2 

97.7 

103.0 

109.2 

113.5 

117.6 

126.2 

129.6 

82 

92.5 

95.3 

98.8 

104.1 

110.4 

114.7 

118.8 

127.4 

130.9 

83 

93.6 

96.3 

99.9 

105.3 

111.6 

115.9 

120.0 

128.6 

132.1 

84 

94,7 

97.4 

101.0 

106.4 

112.7 

117.1 

121.2 

129.9 

133.4 

85 

95.7 

98.5 

102.1 

107.5 

113.9 

118.3 

122.4 

131.1 

134.6 

86 

96.8 

99.6 

103.2 

108.7 

115.0 

119.4 

123.6 

132.3 

135.9 

87 

97.9 

100.7 

104.3 

109.8 

116.2 

120.6 

124.8 

133.6 

137.1 

88 

98.9 

101.7 

105.4 

110.9 

117.4 

121.8 

125.9 

134.8 

138.4 

89 

100.0 

102.8 

106.5 

112.0 

118.5 

123.0 

127.1 

136.0 

139.6 

90 

101.1 

103.9 

107.6 

113.1 

119.7 

124.1 

128.3 

137.3 

140.9 

91 

102.1 

105.0 

108.7 

114.3 

120.8 

125.3 

129.5 

138.5 

142.1 

92 

103.2 

106.1 

109.8 

115.4 

122.0 

126.5 

130.7 

139.7 

143.4 

93 

104.2 

107.1 

110.8 

116.5 

123.1 

127.7 

131.9 

141.0 

144.6 

94 

105.3 

108.2 

111.9 

117.6 

124.3 

128.8 

133.1 

142.2 

145.8 

95 

106.4 

109.3 

113.0 

118.8 

125.4 

130.0 

134.3 

143.4 

147.1 

96 

107.4 

110.4 

114.1 

119.9 

126.6 

131.2 

135.5 

144.6 

148.3 

97 

108.5 

111.4 

115.2 

121.0 

127.7 

132.3 

136.7 

145.9 

149.5 

98 

109.5 

112.5 

116.3 

122.1 

128.9 

133.5 

137.8 

147.1 

150.8 

99 

110.6 

113.6 

117.4 

123.2 

130.0 

134.7 

139.0 

148.3 

152.0 

100 

111.7 

114.7 

118.5 

124.3 

131.2 

135.8 

140.2 

149.5 

153.2 

105 

116.9 

120.0 

123.9 

129.9 

136.9 

141.6 

146.1 

155.6 

159.4 

no 

122.2 

125.4 

129.4 

135.5 

142.6 

147.4 

152.0 

161.6 

165.5 

115 

127.5 

130.7 

134.8 

141.0 

148.3 

153.2 

157.8 

167.7 

171.6 

120 

132.8 

136.1 

140.2 

146.6 

153.9 

159.0 

163.7 

173.7 

177.7 

125 

138.1 

141.4 

145.6 

152.1 

159.6 

164.7 

169.5 

179.7 

183.7 

130 

143.3 

146.7 

151.0 

157.6 

165.2 

170.4 

175.3 

185.6 

189.8 

135 

148.6 

152.0 

156.4 

163.1 

170.9 

176.2 

181.1 

191.6 

195.8 

140 

153.9 

157.4 

161.8 

168.6 

176.5 

181.9 

186.9 

197.5 

201.8 

145 

159.1 

162.7 

167.2 

174.1 

182.1 

187.6 

192.6 

203.4 

207.7 

150 

164.3 

167.9 

172.6 

179.6 

187.7 

193.2 

198.4 

209.3 

213.7 

160 

174.8 

178.6 

183.3 

190.5 

198.9 

204.6 

209.9 

221.1 

225.5 

170 

185.3 

189.1 

194.0 

201.4 

210.0 

215.8 

221.3 

232.8 

237.4 

180 

195.7 

199.7 

204.7 

212.3 

221.1 

227.1 

232.7 

244.4 

249.1 

190 

206.2 

210.2 

215.4 

223.2 

232.2 

238.3 

244.0 

256.0 

260.8 

200 

216.6 

220.7 

226.0 

234.0 

243.2 

249.5 

255.3 

267.6 

272.5 

220 

237.4 

241.8 

247.3 

255.6 

265.2 

271.7 

277.8 

290.6 

295.7 

240 

258.2 

262.7 

268.5 

277.2 

287.1 

293.9 

300.2 

313.5 

318.8 

260 

279.0 

283.7 

289.6 

298.6 

308.9 

316.0 

322.5 

336.2 

341.7 

280 

299.7 

304.5 

310.7 

320.0 

330.7 

338.0 

344.7 

358.9 

364.5 

300 

320.4 

325.4 

331.8 

341.4 

352.4 

359.9 

366.9 

381.5 

387.3 
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TABLE VIII.- DECIBEL SPREAD FOR SPECIFIED 
CONFIDENCE PERCENTAGES 







TSA Output 
Calcomp Plot 
Computer 
Program 


Figure 1.- Dynamic data reduction process. 
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Figure 2.- Simplified flow through TSA program. 







Gain 



(a) Low pass. 

Figure 3.- Precondition simple filters. 
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Normalized Frequency, f 

(b) High pass. 

Figure 3.- Continued. 



-2 


3 


Normalized Frequency, f 

(c) Band pass 
Figure 3.- Continued. 
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(d) Band reject. 
Figure 3.- Continued. 
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Gain 


1 -00 



(e) First difference. 
Figure 3.- Concluded. 
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(a) First difference and band reject cascade. 
Figure 4.- Examples of precondition compound filters. 
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(b) Band reject and band reject cascade. 
Figure 4.- Concluded. 
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Figure 5.- Recursive filter. 
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(c) Arbitrary analysis. 


Figure 6.- Examples of the analysis band -pass filter options. 
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(c) Second term. Third term. 

Figure 9.- Hermite functions. 
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(q) Sixteenth term. 


(r) Seventeenth term. 




Figure 9 .- Concluded. 
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